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Final  Report  for  AFOSR  89-0323 

Work  has  continued  on  analysis  of  models  of  fluid  dynamic  phenomena.  A  number  of  specific 
e.xamples  are  given  in  what  follows.  In  some  of  these  work  is  in  progress. 

Section  1.  Applications  of  free  streamline  modeling. 

a)  Interaction  with  viscous  boundary  layers. 

The  steady  flow  past  a  lifting  surface  at  high  Reynolds  numbers  and  modest  angle  of  attack  can  be 
thought  of  as  consisting  of  two  parts,  the  exterior  flow  in  which  the  flow  is  essentially  inviscid,  and  a 
thin  region  near  the  body  in  which  viscosity  is  important,  the  “boundary  layer.”  In  this  boundary 
layer  special  assumptions  about  the  nature  of  the  flow  and  the  thinness  of  the  domain  allow 
approximations  which  simplify  the  equations  to  be  solved  considerably.  If  the  simplified  equations 
have  been  solved  the  mathematical  problem  becomes  one  of  matching  the  solution  with  an  external 
potential  flow.  These  ideas  go  back  to  Prandtl  [P].  Over  the  years  the  analysis  of  the  boundary  layer 
equations  has  been  refined  by  many  workers.  In  particular  for  many  two  dimensional  problems 
averaged  equations  which  lead  to  a  system  of  ordinary  differential  equations  for  averaged  quantities 
have  been  given.  The  inclusion  of  transition  from  laminar  to  turbulent  flow,  and  the  continuation  of  a 
turbulent  boundary  layer  can  also  be  included  if  certain  empirical  relations  are  accepted.  The  primary 
quantities  solved  for  in  these  equations  are  the  momentum  thickness  and  displ2w:ement  thickness;  these 
quantities  may  be  thought  of  as  describing  the  development  of  the  boundary  layer  downstream  along 
the  body.  (A  description  of  one  version  of  these  equations,  is  be  given  below.  General  references  to 
this  subject  are  [Sch]  and  [M],  chapter  5.) 

The  body  plus  the  boundary  layer  “strip”  may  be  thought  of  as  a  domain  complementing  a 
potential  flow  region.  Since  the  boundary  layer  equations  require  the  pressure  distribution  on  the 
surface  (from  an  exterior  potential  flow)  as  data,  an  iteration  is  suggested  in  which  the  original  body  is 
augmented  by  the  displacement  thickness,  a  new  potential  flow  is  computed,  the  solution  of  the 
boundary  layer  equations  with  the  new  pressure  distribution  is  obtained,  and  so  on.  In  fact,  various 
^versioq  of  this  have  been  in  long  use  in  the  aerodynamics  community. 

If  the  tangential  component  of  velocity  in  the  boundary  layer  becomes  negative,  the  layer  will 
thicken  and  in  a  short  distance  the  assumptions  leading  to  the  boundary  layer  equations  are  no  longer 
satisfied.  One  says  the  boundary  layer  has  “separated.”  If  separation  occurs  a  global  description  of  the 
flow  along  the  lines  of  the  iteration  alluded  to  above  is  much  more  difficult.  In  this  proposal  we  will 
attempt  to  deal  with  this  problem.  We  combine  a  boundary  layer  calculation  with  a  free  streamline 
wake  flow  for  the  outer  potential  flow.  In  the  iteration,  a  separation  point  is  determined  by  the 
boundary  layer  equations,  a  wake  flow  is  computed  and  the  boundary  layer  equations  are  solved  again 
with  the  new  pressure  distribution  on  the  unseparated  part  of  the  boundary.  The  wake  flow  uses 
Tulin’s  double  spiral  vortices  to  terminate  a  near  wake  and  connect  with  free  streamlines  which 
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come  together  at  infinity  downseam.  This  model  was  proposed  and  developed  in  Lpfef];  a  version  with 
improved  numerical  implementation  is  described  below.  The  body  is  approximated  by  a  polygon,  and 
Schwarz-Christoffel  methods  are  used;  the  computational  work  is  in  computing  Schwarz- 

Christoffel  parameters  and  once  this  is  done  the  solution  is  given  by  an  analytic  expressioi^ 

It  was  shown  in  [ET]  by  Elcrat  and  L.N.  Trefethen  that,  in  the  classical  model  of  a  wake  originally 
formulated  by  Helmholtz  [H]  and  KirchhofT  [K],  the  flow  past  a  polygonal  obstacle  can  be  given  in 
terms  of  a  single  analytic  formula.  Further,  the  disposable  parameters  can  be  computed  reliably  and 
efficiently  using  an  adaptation  of  the  earlier  work  of  Trefethen  [Tre]  on  the  Schwarz-Christoffel 
transformation.  If  =  u  —  iv  where  the  velocity  components  are  u,  v,  and  an  upper  half  i-plane  is 
used  as  a  parameter  domain,  the  velocity  in  the  Helmholtz-Kirchhoff  flow  can  be  written 

where 

'■k  =  ((<  - 't)/  (1  -  -  \/i^  )f'‘- 

In  this  formula  =  Tfc  -t- 1  ~  Tfc  where  7^  is  the  angle  the  fcth  side  makes  with  the  horizontal,  is 
the  parameter  value  corres  iding  to  the  vertex  Zj^,  and  •  denotes  the  set  ...,  f*,  where 

corresponds  to  the  stagnation  point  z^. 

The  formula  for  the  flow  is 

dz  _  dw 
dt  ^0 

1  0 

where  tv  =  ^  W(<  —  is  the  complex  velocity  potential,  and  W  is  a  positive  constant. 

The  parameter  is  determined  in  terms  of  fj,  ...,  by  Co(°°)  = 

.  n-1 

<*  =  -  cos  0k  (  “  ^i))- 

A  system  of  equations  for  W,  ...,  is  obtained  by  imposing  the  side  length  conditions;  the  side 

lengths  are  computed  using  compound  Gauss-Jacobi  quadrature  as  in  [ETj.  The  equations  can  be 
solved^in  a  fast  and  efficient  manner  for  non-pathological  geometries  with  a  modest  number  of  sides 
(up  to  30  in  our  experience)  using  a  Newton  type  method  with  the  equally  spaced  values  as  an  initial 
guess.  (In  [Tre]  and  [ET]  the  routine  NSOIA  of  Powell  [PO]  was  used;  in  related  work  Dias  [DET], 
[DI]  used  the  IMSL  routine  ZXSSQ.) 

The  Helmholtz-Kirchhoff  wake  model  suffers  from  the  diseidvantage  that  the  constant  pressure  in 
the  “deadwater”  region  between  the  free  streamlines  must  equal  the  free  stream  pressure;  in 
experimentally  realized  flows  the  pressure  may  be  nearly  constant  but  is  usually  lower  than  the  free 
stream  pressure.  Unfortunately,  it  is  mathematically  impossible  to  have  free  streamlines  that  join 
together  at  a  finite  point  with  a  constant  pressure  unless  it  is  larger  than  the  free  stream  value.  (This 
is  the  Brillouin  paradox;  see  Gilbarg  [G]  for  example.)  This  has  lead  to  a  number  of  attempts  to 
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“close”  the  wake  with  some  artificial  device  in  order  to  have  an  underpressure  and  get  a  reasonable 
description  of  the  flow  near  the  body.  (Several  of  these  are  discussed  in  [G];  see  also  [BZ],  [GU].)  We 
are  interested  here  in  one  of  these  in  particular,  the  double  spiral  vortex  model  of  Tulin  [T].  In  this 
model  there  is  a  jump  discontinuity  in  the  pressure  up  to  the  free  stream  value  at  some  point 
downstream  on  each  free  streamline.  By  Bernoulli’s  equation  this  is  equivalent  to  a  jump  in  speed 
from  a  value  q  to  the  free  stream  value  which  we  take  to  be  1. 

By  necessity  the  argument  of  the  velocity  must  become  infinite  at  such  a  point  and  the  left  and 
right  branches  spiral  around  each  other  infinitely  many  times  in  approaching  each  other,  hence  the 
name  given  to  the  model.  Tulin  used  these  ideas  for  a  linearized  wake  model,  but  “full  nonlinear” 
versions  have  Keen  implemented.  (See  [Fu]  for  example.)  A  version  of  this  was  given  for  polygonal 
obstacles  in  [BEl].  We  introduce  the  parameters  tg  >  I  and  tg/  <  —  1  and  the  functions 


Xs(t)  =  {(t  -  tg)/(ttg  -  1  +  yijn  ))■'" 

and 

Xg,{i)  =  (( 1  -  ttg>  +  /  (f  _ 

Then  ^  Xg  Xgi  will  be  the  complex  velocity.  The  equation  for  is  replaured  by  a  new 

expression  involving  tg,  tgi  also.  In  order  to  have  a  uniquely  determined  flow  two  additional 
conditions  are  needed  to  determine  tg,  tgi  .  In  [BEl]  we  used  two  algebraic  equations  involving  the 
parameters  which  were  derived  from  the  condition  that  the  asymptotic  width  of  the  far  wake  was  zero. 
A  number  of  successful  computations  were  done  in  the  work  leading  to  [BEl],  but  the  side  length 
conditions  coupled  with  these  algebraic  equations  are  a  poorly  conditioned  system  and  the 
computational  effort  was  much  greater  than  that  in  [ET].  This  has  led  us  to  reformulate  the 
computational  problem  as  follows.  Suppose  that  Zg  and  Zgi  are  the  spiral  points.  Then  we  require 
that 

R€al{zg-Zg,)  =  0, 

(the  near  wake  has  a  symmetric  closure),  and  that 

Im{z(T)-z{-T))  =  Q 

for  some  “large”  value  of  T.  This  second  condition  is  a  finite  approximation  to  the  condition  of  zero 
asymptotic  width  for  the  far  wake.  (In  an  airfoil  problem  a  value  of  J’=  4  typically  puts  z{T)  about 
20  chord  lengths  downstream.)  These  are  then  two  additional  conditions  and  the  count  is  right  in 
determining  W,t^,  ...,  tg,  tgi.  The  new  equations  involve  integrals  which  are  evaluated  in 

essentially  the  same  way  as  the  others,  and  the  system  of  equations  obtained  is  well  behaved.  When 
both  methods  work  (the  present  method  and  the  one  presented  in  [BEl]),  the  results  differ  by  less  than 
a  few  percent. 

We  give  here  the  boundary  layer  equations  in  the  form  that  we  will  use  then. 

The  motion  of  a  viscous  fluid  flowing  in  a  thin  layer  near  a  solid  surface  can  be  described  by  the 


equations 


uux  +  vuy=  -  Px  +  ^yy 

Ux  +  Vy  =  0 

where  u  =  i;  =  0ony  =  0,  (u,  v)  approaches  (t/ (l),  0)  as  y  goes  to  infinity,  and 

p(x)  +  i  U‘^{x)  =  constant. 

The  equations  have  been  written  in  dimensionless  form  and  R  is  the  Reynolds  number;  curvature  of 
the  surface  has  been  neglected.  The  solutions  of  these  equations  are  an  asymptotic  limit  of  solutions  of 
the  Navier-Stokes  equations  which  are  to  be  matched  with  an  outer  potential  flow  which  has  velocity 
distribution  U (x)  on  the  edge  of  the  boundary  layer.  In  particular,  the  solutions  may  be  desired 
downstream  of  a  stagnation  pioint.  Then  a  useful  approximation  can  be  obtained  by  averaging  to 
obtain  ordinary  differential  equations  for  the  displacement  thickness, 

00 

5(x)  =  J (I -^)dy 

and  momentum  thickness,  0 

oo 

^(x)=  J  ^{l-^)dy. 

0 

If  the  flow  is  turbulent  one  set  of  such  equations  are 

6' +  {H+2)  {U'/U)e  =  \c^, 

(the  momentum  integral  equation),  and 

^  =  .0306  (Hj 

where  H=6/B  is  the  “shape  factor”  and  is  the  turbulent  -^kin  friction.  We  will  take  to  be  given 
by  the  Ludwieg-Tillman  expression 

Cf  =  .246  X  (iJ(7^)--268 

is  an  (empirically  derived)  function  of  H{c{.  Moran  [M]). 

In  the  initial,  laminar,  flow  region  we  use  the  method  of  Thwaites.  Here  an  empirical  correlation 
between  H  and  the  skin  friction  allows  exact  integration  of  the  momentum  integral  equation. 
Transition  is  determined  by  Michel’s  method,  and  after  transition  the  above  equations  (Head’s  method) 
are  used.  Our  implementation  is  essentially  the  program  INTGRL  given  by  Moran  [M].  This  code  is 
used  as  a  subroutine  which  integrates  downstream  from  a  stagnation  point  until  separation  occurs 
(either  laminar  or  turbulent).  In  the  turbulent  region  separation  is  taken  to  have  occurred  when  H 
reaches  a  certain  value.  In  INTGRL  this  value  is  2.4;  we  have  used  the  value  1.8.  The  value  we  have 
taken  was  observed  in  one  set  of  experiments  for  a  GA(W)-2  airfoil  that  we  discuss  below,  and  this 
choice  led  to  excellent  agreement  with  these  experiments. 

We  will  describe  now  in  more  detail  the  iteration  that  we  have  in  mind  for  separated  flows. 
Suppose  that  at  a  given  angle  of  attack  a  potential  flow  (without  separation)  has  been  obtained  for  the 


obstacle  under  consideration,  and  that  the  pressure  distribution  is  known  as  a  function  of  arc  length 
downstream  from  the  stagnation  (loint.  (We  are  tacitly  ctssuming  an  airfoil  section  at  positive  angle  of 
attack  so  that  only  separation  on  the  top  is  investigated.  This  is  not  essential  and  the  generalization  of 
this  procedure  to  other  situations  is  straightforward).  This  is  then  used  as  data,  together  with  the 
Reynolds  number,  for  a  boundary  layer  calculation,  the  output  of  which  is  the  separation  point,  if  any, 
and  9  and  H  as  functions  of  arc  length  up  to  separation.  If  separation  occurs  a  wake  flow  is  computed 
as  described  in  section  1  using  the  separation  point  and  separation  pressure  from  the  boundary  layer 
calculation  as  data.  The  displacement  thickness  is  added  to  the  current  profile.  This  provides  a  new 
pressure  distribution  which  can  be  used  as  input  to  the  boundary  layer  subroutine.  The  process  is 
iterated  until  the  separation  point  does  not  move  further  upstream. 

We  report  here  the  results  of  some  numerical  calculations  and  make  comparisons  with 
experiments. 

In  these  examples  we  have  used  as  input  the  pressure  distribution  for  attached  flows.  Typically 
convergence  was  obtained  for  a  given  angle  of  attack  in  2-3  iterations  with  a  total  CPU  time  of  less 
than  2  minutes  on  the  IBM  3  ES-9000,  model  440  at  Wichita'State  University. 

The  experiments  of  Wentz,  Seetheram,  and  Rogers  reported  in  [WSR]  were  our  principal 
benchmark.  These  were  carried  out  at  a  Reynolds  number  of  2.2  million,  and  we  show  in  figures  1-3 
below  our  computations  plotted  with  experimental  values  for  several  angles  of  attack  for  a  GAW-2 
airfoil. 


In  addition  we  have  run  the  case  a=18”,  R—  4.3  x  10®  and  compared  with  the  experimental  curve  given 
in  [MB].  The  result  is  shown  below,  in  figure  4.. 


rituva  3. 

A  major  impetus  for  this  work  is  its  use  in  matching  for  3-D  flows  as  in  [BE2].  Work  is  in 
progress  on  doing  this  when  separation  on  a  section  is  determined  by  the  Reynolds  number.  The 
matching  leads  to  a  fixed  point  equation  relating  the  parameters  in  the  inner  flow  on  a  two 
dimensional  cross  section  to  the  circulation  and  source  strength  in  a  lifting- line-source  model  of  the 
outer  solution.  The  velocity  potential  of  the  outer  flow  has  an  inner  expansion 

4>  =  x-y-^  j  ^  /l(2)  arctan  {y/x)  -f-  i  f2{z)  in  r -h  0(e) 

where  z  is  transverse  coordinate,  £  is  the  aspect  ration,  and  the  integral  is  a  Cauchy  principal  value. 
The  inner  solutions  has  an  outer  expansion 

$  =  I  -  a(z)  y  -f  B{z)  Inr  -  C{z)  arctan  [y/x)  -f  0(e) 
where  a  is  the  downwash.  Matching  to  first  order  in  £  requires 
fliz)  =  2xC{z),  /2(^)  =2  3rB(2). 

The  downwash  is  then  given  by 


B{z),  C{z)  are  functions  of  the  parameters  in  the  wake  flow,  and  this  yields  a  fixed  p>oint  equation 
which  can  be  solved  numerically  as  in  (BE2]. 

The  separation  point  and  pressure  on  a  section  are  determined  cis  described  previously,  so  that  the 
model  obtained  of  the  3-D  wake  contains  cross  flow  pressure  gradients  and  is  not  a  constant  pressure 
region  as  in  [BE2].  We  can  formulate  a  crude  model  of  the  flow  in  the  3-D  wake  with  changes  in  the 
sign  and  magnitude  of  B{z)  (this  is  essentially  the  source  strength);  there  is  some  hope  that  we  may  be 
able  to  use  this  to  help  explain  the  recirculation  cells  in  the  wake  behind  a  wide,  bluff  body.  Work  is 
in  progress  on  this. 

b)  Mixed  design-analysis  for  2-D  airfoils. 

In  the  work  described  above  for  separated  flow  on  an  airfoil  it  was  found  that  a  polygonal 
approximation  to  an  airfoil  provides  surprisingly  good  approximations  for  pressure  distributions.  This 
requires  some  explanation. 

The  exact  velocity  profile  (and  hence  pressure  distribution)  of  a  flow  past  a  polygon  has  infinities 
(or  zeros)  at  each  vertex.  It  is  extremely  ragged,  nonmonotone  and  a  poor  approximation  to  that  of  a 
smooth  body.  On  the  other  hand  if  pressures  are  taken  at  mid-pK>ints  of  straight  boundary  segments 
and  smoothly  interpolated  this  difficulty  disappears.  In  fact  this  is  closely  related  to  common  practice 
with  panel  methods  in  aerodynamics.  Our  recent  experience  has  been  that  a  Schwarz-Christoffel 
solution  does  surprisingly  well  with  a  relatively  small  number  of  sides  (“panels”). 

Since  a  wake  model  is  nothing  but  a  replacement  of  flow  tangency  with  a  specification  of  the 
pressure  over  a  part  of  the  flow,  the  suggestion  is  natural  that  this  can  be  done  in  a  more  general  way 
than  simply  requiring  the  pressure  to  be  constant  on  the  wake  streamlines.  To  sjjecify  pressures  over  a 
part  of  the  airfoil  or  wake  boundary  leads  to  a  “generalized  parameter  problem”  for  a  Schwarz- 
Christoffel  transformation,  as  described  in  [Tr  RD].  An  interesting  flow  problem  which  can  be  used  as 
a  first  pass  in  testing  the  feasibility  of  this  idea  is  the  formation  of  a  “leading  edge  separation  bubble” 
as  depicted  below. 


We  will  begin  by  investigating  this  problem, 
c)  Crack  detection. 

Another  problem  that  we  are  working  on  is  described  below.  This  is  not  a  fluids  problem.  It  is 
placed  here  in  the  discussion  because  of  the  methods  to  be  used. 


Consider  a  polygon  P  with  four  distinguished  vertices  a,b,c,d.  If  fixed  constant  voltages  are 
imposed  on  the  arcs  ab  and  cd  and  the  rest  of  the  boundary  is  insulated,  the  resistance  of  this 
polygonal  “circuit  element”  can  easily  be  measured  experimentally  or  computed  using  the  solution  of  a 
mixed  boundary  vlaue  problem  for  Laplace’s  equation.  The  computation  has  been  give  a  particularly 
elegant  formulation  by  Trefethen  in  [TrRD]  using  the  Schwarz-Christoffel  transformation.  Suppose 
now  that  a  crack  of  unknown  (p>olygonal)  shape  extends  into  P.  If  this  crack  is  described  by  N 
geometrical  parameters  then  an  algorithm  for  their  determination  can  be  given  as  follows.  Suppose 
that  the  “terminals”  ab  and  cd  are  moved  around  the  boundary  N  times  to  give  N  different  values 
for  the  resistance.  These  measurements  then  give  N  conditions  with  which  to  determine  the  geometry 
of  the  crcick.  Trefethen  did  essentially  this  in  the  case  iV=l  although  he  considered  the  crack  length  as 
a  variable  to  be  tuned  in  order  to  obtain  a  prescribed  resistance.  The  determination  of  several 
geometrical  parameters  is  another  “generalized  parameter  problem”  for  a  Schwarz-Christoffel  map  and 
it  seems  very  likely  that  this  proposed  algorithm  can  be  made  to  work. 

A  similar  idea  can  be  put  forward  if  the  crack  is  interior  to  the  domain.  The  physics  is  the  same, 
but  a  doubly  connected  version  of  the  polygonal  mapping  problem  must  be  used.  Fortunately,  the 
computational  details  of  this  have  been  worked  out  by  Daeppen  [Daep].  The  program  that  he  wrote  is 
not  available,  but  we  plan  to  write  our  own  program  based  on  his  ideas  as  a  part  of  this  project.  For  a 
three  dimensional  problem  conformal  mapping  is  not  available  and  a  computational  algorithm  will  be 
more  difficult,  but  is  seems  likely  that  the  two  dimensional  problem  will  serve  as  a  guide  for  what  is 
possible. 

Section  2.  Variational  Methods  tor  Vortex  Flows. 

a)  Steady  Flows. 

The  use  of  variational  methods  for  the  study  of  steady,  planar  flows  of  an  ideal  fluid  with  regions 
of  nonzero  vorticity  has  been  in  progress  since  the  time  of  Kelvin  [K].  More  recently,  Benjamin  [B76], 
Turkington  and  Friedman  [FT],  and  Burton  [BVR]  have  studied  vortex  rings  using  this  approach.  The 
planar  problem  was  studied  by  Turkington  [T83],  Burton  [Bj,  and  Elcrat  and  Miller  [EMI],  [EM2]. 
Numerical  results  have  been  given  in  [EM3]  and  [TE].  The  work  of  [B]  is  of  particular  importance  for 
this  proposal  in  that  the  vorticity  in  the  ”core”  regions  need  not  be  constant;  the  class  in  which 
variational  solutions  are  found  is  the  set  of  rearrangements  of  a  fixed  function  (ie.  functions  whose  level 
sets  have  the  same  measure)  and  the  profile  function  which  relates  the  vorticity  to  the  stream  function 
is  determined  as  a  part  of  the  solution.  (This  is  in  line  with  Benjamin’s  original  suggestions  in  [B76].) 
One  difficulty  with  Burton’s  theory  is  that  it  does  not  always  apply  directly  to  situations  of  physical 
interest.  We  plan  to  extend  his  results  to  flow  in  a  bounded  domain  with  inflow  and  outflow  boundary 
conditions  on  parts  of  the  boundary.  This  can  be  used  to  model  flow  in  a  channel  as  well  as 


approximating  the  flow  past  an  obstacle  in  an  unbounded  domain.  As  a  part  of  our  efforts  we  plan  to 
develop  an  iterative  scheme  for  the  solution  in  which  the  profile  function  and  hence  the  vorticity  are 
determined  iteratively  using  variational  methods.  After  successful  completion  of  this  task,  a  r  umerical 
solution  of  the  problem  can  be  obtained  in  a  straightforward  way.  In  a  related  work  axisymmetric 
flows  past  symmetric  bodies  will  be  studied. 

These  flow  problems  are  relevant  to  high  Reynolds  number  flows  past  bluff  bodies  and  in  the 
interior  regions  such  as  an  air  duct.  In  practical  situations  the  flow  in  these  configurations  is  likely  to 
be  unsteewiy  and  turbulent,  but  steady  flows  can  be  useful  in  approximating  averaged  quantities;  also 
the  possibilty  of  introducing  external  controls  which  maintain  the  stability  of  a  ste£idy  flow  is  a  Holy 
Grail  for  fluid  dynamicists,  and  further  information  about  steady  flows  improves  the  chances  for 
attaining  this  goal. 

The  iterative  scheme  mentioned  above  will  now  be  explained  in  the  context  of  the  model  problems 
of  maximizing  the  energy  for  vortices  whose  support  are  contained  in  a  fixed  bounded  domain.  A 
maximum  is  sought  for 

E(u) 

where  K  is  the  Green’s  operator  for  the  negative  Laplacian  in  fi,  and  the  maximum  is  sought  in  the 
class  of  functions  which  are  rearrangements  (level  sets  have  the  same  measure)  of  a  fixed  function  in 
L^(Cl).  (  If  the  stream  function  rp=Ku  is  introduced,  an  integration  by  parts  shows  that  E(ijJ)  is  the 
Dirichlet  integral  of  xp).  The  iteration  that  we  will  work  with  is  based  on  ideas  given  by  Eydeland  and 
Turkington  [EyTu].  A  sequence  converging  to  a  maximizer  is  generated  by 

uj^  argmax^  J 
Q 

VVe  will  take  the  maximum  over  the  class  of  rearrangments  defining  the  original  variational  problem. 
We  have  proven  that  this  sequence  converges  in  L^{Q,)  and,  further,  that  there  are  increasing  profile 
functions  4>y  such  that  {xp^)  =  w-',  and  that  the  •sequence  (pj  converges  almost 

everywhere.  We  plan  to  give  a  numerical  implementation  of  this  iteration  by  approximating  the 
rearrangement  condition  at  a  finite  number  of  values. 

There  is  another  class  of  variational  problems  for  2-D  vortex  flows  that  we  have  considered.  The 
physics  of  these  flows  can  be  associated  with  large  Reynold’s  number  limits  for  steady  solutions  of  the 
Navier  Stokes  equations  in  the  flow  past  a  bluff  body.  (See,  for  example,  [S],  in  particular  the 
discussion  of  “Prandtl-Bachelor”  flows.)  They  can  also  be  thought  of  as  having  practical  value  as 
averages  for  turbulent  flows.  The  problem  is  to  find  a  flow  which  contains  a  vorticial  region  in  the 
separation  zone  behind  an  obstacle,  for  example  as  in 


The  •  <  iational  problem  for  such  flows  has  some  formal  similarity  to  that  discussed  above,  but  is 
different  in  subtle  ways.  In  particular,  if  the  vorticity  is  constant,  say  equal  to  -u;,  then  a  functional 
which  is  associated  with  this  problem  is  given  by 

1  =  5  J  +  u)  J  Tp 

D  B 

where  0  is  the  stream  function,  D  is  a  “control  volume”  in  which  the  flow  takes  place,  and  B  is  the 
vorticial  region.  In  fact,  a  number  of  interesting  flows  were  computed  by  an  interation  which 
minimizes  this  functional  by  M.A.  Goldshtik  several  years  ago  [Gold].  (This  work  does  not  appear  to 
be  widely  known  in  the  west.  VVe  learned  about  it  during  a  recent  visit  by  Goldshtik  to  Wichita.)  A 
surprising  fact  about  this  functional  is  that  a  maximizer  B,  if  it  existed,  could  not  have  interior  points. 
The  functional  that  was  considered  by  Turkington  in  [T83],  for  example,  can  be  brought  into  the  above 
form:  in  that  work  a  maximizer  was  found  subject  to  a  constraint  on  the  area  (circulation)  of  the 
vortex  patch. 

We  plan  to  pursue  further  work  along  these  lines  using  the  ideas  of  Goldshtik.  In  particular  a 
“vortex  sheet”  or  jump  in  the  speed  across  the  boundary  of  the  vortex  region  arises  if  the  term  0'\B\  is 
subtracted  from  /;  <7  is  the  magnitude  of  the  jump.  We  have  given  a  theoretical  treatment  for  a  model 
problem  in  which  a  vorticial  region  bounded  by  a  vortex  sheet  is  contained  in  a  bounded  domain.  The 
minimizer  is  bounded  by  a  Lipschitz  curve,  and  an  essential  part  of  the  analysis  is  to  establish 
Hadamard  typw  variational  formulas  (Garabedian,  PDE,  chatper  15)  for  Lipschitz  domains.  We  hope 
to  give  a  similar  treatment  to  problems  of  the  type  considered  by  Goldshtik. 

b)^  Hamiltonian  formulation  of  time  dependent  wakes. 

Consider  a  two  dimensional  flow  in  which  there  are  bounded  regions  in  which  the  fluid  motion  has 
constant  vorticity.  The  flow  can  be  described  by  harmonic  velocity  potentials  in  the  various  regions; 
these  are  perturbations  from  a  fixed  constant  vorticity  flow  in  the  rotational  regions.  A  complete 
description  of  the  flow  requires  these  time  dependent  potentials  and  the  boundaries  of  the  rotational 
regions. 

To  fix  ideas  consider  a  symmetric  flow  past  a  symmetric  body  which  is  uniform  at  infinity  and  has 
a  single  vortex  region  in  the  upper  half  plane  adjacent  to  the  body.  The  situation  is  depicted  below. 


n 


In  Q  the  velocity  is  given  by  SP  (p  +  Ui  and  in  by  V  0*  +  yC  —  y,  x)  +  Ui={ipy,  —  ipx)- 
In  recent  work  we  have  shown  that  this  problem  can  be  given  a  Hamiltonian  formulation.  This  uses 
ideas  recently  develop>ed  by  Benjamin  [BeBr],  [Be],  [BWWj.  If  K  is  the  perturbation  energy  in  the 
flow  then 

A'  =  1  V  ci|2  +  1  y  (I  V  ,^'|2  +  _  3/,  .  V  <&'  +  /) 

"n  Q' 

-U 

J 

where  f  is  a.  fixed  function,  we  have  taken  the  density  of  the  fluid  to  be  1,  and  $=  (f>(x,h{x,t),  f) 
h(x,  t),  t).  We  denote  by  SK/Su  a  variational  derivative  with  respect  to  a  quantity  u  as  in 
[BeBr],  [BWW]. 

Let  4*  =  ^  —  We  have  shown  then  that 
A  -6K  r  - 

These  express  kinematical  compatibility  and  the  Bernoulli  equations  in  H,  We  plan  to  extend  this 
to  nonsymmetric  flows,  multiple  vorte.x  configurations,  and  axisymmetric  flows. 

There  are  several  important  features  of  this  new  formulation  of  the  unsteady  wake  problem.  The 
evolution  is  described  in  terms  of  quantities  defined  only  on  the  boundary  of  the  vorticial  region 
bounded  by  vortex  sheets,  but  the  formulation  is  not  inherently  singular  as  is  the  case,  for  example, 
with  the  BirkofT-Rott  equation  for  evolution  of  vortex  sheets;  a  desingularization  should  not  be 
necessary  for  a  numerical  solution  as  in  [Kr].  Accurate  and  efficient  computation  of  Dirichlet  Integrals 
of  harmonic  functions  is  required,  of  course,  but  we  have  extensive  experience  in  doing  this,  and  work  is 
in  progress  in  further  developing  this  capability  (cf.  discussion  of  conformal  mapping  below.) 

The  Hamiltonian  formulation  leads  naturally  to  variational  principles  for  steady  flows.  In 
particular,  the  computation  of  Prandtl- Bachelor  flows  and  the  establishment  of  their  relevance  for  large 
Reynolds  number  viscous  flows  is  a  long  standing  problem  ([S]),  and  the  new  variational  principles  that 
arise  from  this  work  may  be  of  help  in  resolving  this. 

If  the  wake  problem  is  idealized  by  removing  the  obstacle  (eg.  it  is  far  upstream  and  its  effects  are 
contained  in  the  vorticity  distribution),  a  flow  problem  which  has  much  in  common  with  the  water 
wave  problem  is  obtained  (cf.[BVVW].)  Then  not  only  steady  solutions,  but  steadily  translating 


configurations  can  be  studied.  This  is  relevant  to  unsteady  oscillations  in  the  wake  behind  a  bluff 
body. 

c)  Conformal  Mapping 

The  conformal  map  of  a  canonical  domain  onto  one  with  an  irregular  boundary  is  a  very  useful 
tool  in  solving  partial  differential  equation  problems  involving  harmonic  functions.  If  the  canonical 
domain  is  the  unit  disk  and  the  curve  bounding  the  irregular  domain  is  smooth  the  methods  which  use 
Fourier  series  and  conjugation  to  determine  the  boundary  correspondence  function,  and  hence  the  map, 
are  the  most  widely  used  and  successful.  Among  these  two  that  compute  the  boundary  correspiondence 
s=s(0),  S=arc  length  along  the  target  curve,  using  Newton’s  method  for  updating  the  approximation 
to  s(5)  were  developed  by  Fornberg  [FO]  and  Wegmann[W].  They  may  be  thought  of  as  arising  from 
a  singular  integral  equation  for  s{9)  on  the  unit  circle,  and  the  linear  equation  which  has  to  be  solved 
at  each  Newton’s  step  requires  that 
(4.1)  2(s„(^))  +  AsM  Un(0) 

extend  to  an  analytic  function  inside  the  disk.  Here  2(s)  is  the  coordinate  function  for  the  curve  and 

^  +  Un-  There  are  some  differences  in  detail  in  how  this  problem  is  solved  at  the  discrete 

level  in  the  two  approaches;  Wegmann  solves  a  Riemann- Hilbert  problem  in  the  disk,  whereas  Fornberg 
projects  onto  functions  analytic  in  the  disk.  Both  methods  are  very  effective;  the  nonlinear  iterations 
usually  converge  from  an  arbitrary  initial  guess  although  known  theoretical  results  show  only  local 
convergence.  The  linear  equations  are  solved  using  FFT’s  and  the  solution  is  fast.  (We  have  found 
Fornberg  to  be  more  robust;  Wegmann  requires  fewer  FFT’s  per  Newton  iteration.)  The  method 
provides  as  part  of  the  solution  an  approximate  Taylor  series  for  the  map. 

A  typical  application  was  dealt  with  by  Fornberg,  namely  the  time  evolution  of  a  free  surface 
under  the  influence  of  gravity.  A  new  conformal  map  is  computed  at  each  time  step.  In  this  and  other 
problems  domains  which  are  thin  or  “fingered”  on  part  of  their  boundary  often  arise.  When  this 
happens  the  “crowding”  phenomena  occurs  [G],  [MZ]  and  a  large  number  of  Fourier  coefficients  are 
required.  In  this  case  it  is  natural  to  ask  if  use  of  a  different  canonical  domain  will  make  the  problem 
more  tractible.  In  the  case  of  polygon  a  Schwarz-Christoffel  method  has  been  given  in  [HT].  Another 
new  method  for  doing  this  is  proposed  here. 

The  difficulty  with  Fourier  series  for  thin  domains  is  that  a  large  number  of  terms  in  a  Taylor 
series  are  required  for  resolution,  i.e.  the  problem  is  one  of  discretization  error.  We  will  consider 
domains  which  are  “thin”  in  one  direction  and  use  the  inside  of  an  ellipse  which  is  the  inverse  image 
under  the  Joukowski  map  2=i(f-(~)  of  |<|  =  p  as  our  canonical  domain.  In  place  of  Taylor  series  we 
use  the  Tchebychev  series 

(TS)  5  >^0  +  £  AnTni^y, 

Tfi{z)  are  the  “Faber  polynomials”  for  this  domaifT  and  are  appropriate  for  the  expansion  of  analytic 
functions  here.  It  can  be  shown  that  a  function  defined  on  the  ellipse  decomposes  into  the  sum  of  a 


Faber  series  and  a  function  analytic  on  the  exterior.  The  analogue  of  Fornberg’s  projection,  in  which 
negative  index  Fourier  coefficients  are  set  to  zero,  is  therefore  possible.  It  can  be  expressed  in  the  form: 
integral  of  f  times  equals  zero  for  all  n.  It  is  convenient  to  use  the  Joukowski  transformation  and 
a  scale  change  to  transform  this  to  the  unit  circle.  If  ^  is  the  new  variable  and  we  denote  the 
transplanted  function  again  by  /  this  condition  is  equivalent  to  a.„  =  a„  for  the  Laurent 

coefficients  of  f(Q.  Since  our  ellipse  is  conveniently  parametrized  by  z  =  ^  +  p'^  e'*^)  we  may 

21S  well  think  of  the  equation  (2.1)  as  being  imposed  on  27r  periodic  functions  or  functions  defined  on 
iO 

^=6**'  as  convenience  dictates.  If  P  ^  ,  P  —  are  the  projections  onto  functions  analytic  in  the  interior 
and  exterior  of  the  circle  our  condition  on  the  Laurent  coefficients  takes  the  form 

P_/=IP+  f  +  Jf 

where  J f=  Cq  ,  and  L  multiplies  Taylor  coefficients  by  and  performs  an  inversion.  If  we  write 
the  right  hand  side  of  (4.1)  as  2  +  tU,  we  obtain 

(4.2)  (P_  -IP  +  )  f(7=  -(P_  -IP4.)  z. 

If  this  is  discretized  and  use  is  made  of  the  iV  —  dimensional  discrete  Fourier  transform  we  obtain 

(4.3)  A.  F T  U  —  —  A  Pjy  z 

where  T  is  a.  diagonal  matrix  (v2ilues  of  t)  and  A  has  the  block  structure 


A=  ®  “ 

P  I 


the  (square)  submatrices  having  order  iV/2  with 


o  ... 


The  first  row  arises  from  discretization  of  the  series  (TS)  vanishing  at  the  origin.  This  amounts  to 
iV/2  complex  equations  for  the  N  real  components  of  U,  in  complete  analogy  to  equation  (126)  of 
Fornberg  [FOJ.  We  introduce  the  even-odd  decomposition  as  done  there  and  use  the  Gauss 
decomposition. 


PAr  = 


M 

^l^M 

M 

where  A/=  iV/2  and  W  j  are  diagonal  matrices  with  roots  of  unity  on  the  diagonal.  This  leads  to  an 


equation  of  order  M : 

{P  +  /)  Fr,  To  Uq  +  +  IW^)  FmUi  =  e 

where  e  is  obtained  by  mulitpiying  the  matrices  on  the  left  times  z,  and  Uq,  are  the  even  and  odd 
components  of  U.  Finally,  the  fact  that  U  is  real  must  be  used  to  obtain  a  full  system.  As  in  [FO]  an 
M  X  M  system  can  be  obtained  for  U q,  and  U  ^  \s  obtained  by  back  substitution.  The  matrices 
involved  are  either  unitary  or  have  an  explicit  structure  that  requires  only  0(M)  operations  for 
equation  solving;  a  conjugate  gradiant  step  can  be  done  in  a  comparable  number  of  operations  to  those 
in  [FO).  The  condition  that  a  point  is  fixed  on  the  boundary  of  the  target  curve  can  be  used  as  in  [FO] 
to  compensate  for  the  fact  that  P  + 1  has  an  eigenvalue  near  zero.  Finally,  as  /?— *oo  Fornberg’s 
method  is  retreived,  and  the  eigenvalue  structure  of  the  operators  in  (4.2)  and  matrices  in  (4.3)  can  be 
analyzed  as  in  [WE],  [WI]. 

A  strong  motivation  for  obtaining  this  method  is  to  deal  with  cases  in  which  the  vorticial 
domains  in  inviscid  flow  problems  become  elongated.  In  particular  an  efficient  computation  of  the 
Dirichlet  integral  can  be  done  using  the  transplantation  to  the  ellipse.  This  may  be  of  great  value  in 
both  steady  and  unsteady  problems  in  which  the  Hamiltonian  functional  involves  the  Dirichlet  integral. 
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